/*
 * geometry.cpp
 *
 *  Created on: 2012-12-7
 *      Author: shane
 */

#include "geometry.h"
#include <algorithm>
using namespace std;
namespace EME
{

bool geo_line2::Intersect(const geo_line2 L1, const geo_line2 L2)
{
	return ((max(L1.begin.x, L1.end.x) >= min(L2.begin.x, L2.end.x))
			&& (max(L2.begin.x, L2.end.x) >= min(L1.begin.x, L1.end.x))
			&& (max(L1.begin.y, L1.end.y) >= min(L2.begin.y, L2.end.y))
			&& (max(L2.begin.y, L2.end.y) >= min(L1.begin.y, L1.end.y))
			&& (geo_point2::Multiply(L2.begin, L2.end, L1.begin)
					* geo_point2::Multiply(L2.begin, L2.end, L1.end) <= 0.)
			&& (geo_point2::Multiply(L1.begin, L1.end, L2.begin)
					* geo_point2::Multiply(L1.begin, L1.end, L2.end) <= 0.));
}
bool geo_line2::containPoint(const geo_point2& p)
{
	return fabs(geo_point2::Multiply(begin, end, p)) < distance() * 1e-6
			&& (p.x - begin.x) * (p.x - end.x) <= 0
			&& (p.y - begin.y) * (p.y - end.y) <= 0;
}
geo_object::geo_object() :
		xmin(0.), xmax(0.), ymin(0.), ymax(0.), area(0.)
{
}
geo_object::geo_object(const geo_object& b) :
		xmin(b.xmin), xmax(b.xmax), ymin(b.ymin), ymax(b.ymax), area(
				b.area)
{
}
geo_rectangle::geo_rectangle(const geo_point2& a, const geo_point2&b):geo_object()
{
	if (a.x == b.x || a.y == b.y)
		throw(EME_geometry_logic_error("Rectangle configuration error", 1));
	if (a.x < b.x)
	{
		xmin = a.x;
		xmax = b.x;
	}
	else
	{
		xmin = b.x;
		xmax = a.x;
	}
	if (a.y < b.y)
	{
		ymin = a.y;
		ymax = b.y;
	}
	else
	{
		ymin = b.y;
		ymax = a.y;
	}
	area = (xmax - xmin) * (ymax - ymin);
}
geo_type geo_rectangle::get_type(){return rectangle;}
geo_rectangle::geo_rectangle(const geo_rectangle& b) :
		geo_object(b)
{
}
int geo_rectangle::inside(const geo_point2& p){
	if(p.x>xmax||p.x<xmin||p.y>ymax||p.y<ymin)
		return 0;
	if(p.x<xmax&&p.x>xmin&&p.y<ymax&&p.y>ymin)
		return 1;
	return 2;
}
bool geo_polygon::issimple()
{
	geo_line2 l1, l2;
	for (vector<geo_point2>::size_type i = 0; i < vertexs.size(); ++i)
	{
		l1.begin = vertexs[i];
		l1.end = vertexs[(i + 1) % vertexs.size()];
		int cn = (i == 0 ? vertexs.size() - 3 : vertexs.size() - 3 - i + 1);
		int i2 = 0;
		while (cn)
		{
			l2.begin = vertexs[(i + 2 + i2) % vertexs.size()];
			l2.end = vertexs[(i + 3 + i2) % vertexs.size()];
			if (geo_line2::Intersect(l1, l2))
				break;
			cn--;
			i2++;
		}
		if (cn)
			return true;
	}
	return false;
}
double geo_polygon::area_of_polygon()
{
	if (vertexs.size() < 3)
		throw(EME_geometry_logic_error("Polygon points number < 3", 2));
	double s;
	s = vertexs[0].y * (vertexs.back().x - vertexs[1].x);
	for (vector<geo_point2>::size_type i = 1; i < vertexs.size(); ++i)
		s += vertexs[i].y
				* (vertexs[(i - 1)].x - vertexs[(i + 1) % vertexs.size()].x);
	if (s == 0.)
		throw(EME_geometry_logic_error("Polygon area size is zero", 3));
	return s / 2;
}
bool geo_polygon::isccwise()
{
	vector<geo_point2>::size_type index = 0;
	for (vector<geo_point2>::size_type i = 1; i < vertexs.size(); ++i)
		if (vertexs[i].y < vertexs[index].y
				|| (vertexs[i].y == vertexs[index].y
						&& vertexs[i].x < vertexs[index].x))
			index = i;
	geo_point2 a = vertexs[(index - 1 + vertexs.size()) % vertexs.size()];
	geo_point2 b = vertexs[(index + 1) % vertexs.size()];
	return geo_point2::Multiply(b, a, vertexs[index]) > 0;
}
geo_type geo_polygon::get_type(){return polygon;}
int geo_polygon::inside(const geo_point2& p)
{
	geo_line2 line;
	int n = vertexs.size();
	int count = 0;
	line.begin = p;
	line.end.y = p.y;
	line.end.x = xmin - 1e10;
	for (int i = 0; i < n; ++i)
	{
		geo_line2 side(vertexs[i], vertexs[(i + 1) % n]);
		if (side.containPoint(p))
			return 2;
		if (fabs(side.begin.y - side.end.y) < side.distance() * 1e-6)
			continue;
		if (line.containPoint(side.begin))
		{
			if (side.begin.y > side.end.y)
				++count;
		}
		else if (line.containPoint(side.end))
		{
			if (side.end.y > side.begin.y)
				++count;
		}
		else if (geo_line2::Intersect(line, side))
			++count;
	}
	return count % 2;
}
geo_polygon::geo_polygon(const std::vector<geo_point2>& vertexs_in):geo_object()
{
	vertexs = vertexs_in;
	if (!isccwise())
	{
		reverse(vertexs.begin(), vertexs.end());
	}
	if (!issimple())
		throw(EME_geometry_logic_error(
				"Polygon configuration error, is not simple", 4));
	area = area_of_polygon();
	if (area == 0.)
		throw(EME_geometry_logic_error(
				"Polygon configuration error, area size 0", 5));
	xmin = vertexs[0].x;
	xmax = vertexs[0].x;
	ymin = vertexs[0].y;
	ymax = vertexs[0].y;
	for (vector<geo_point2>::iterator i = vertexs.begin(); i != vertexs.end();
			++i)
	{
		if ((*i).x > xmax)
			xmax = (*i).x;
		if ((*i).x < xmin)
			xmin = (*i).x;
		if ((*i).y > ymax)
			ymax = (*i).y;
		if ((*i).y < ymin)
			ymin = (*i).y;
	}
}
geo_polygon::geo_polygon(const geo_polygon& b) :
		geo_object(b), vertexs(b.vertexs)
{
}
geo_ellipse::geo_ellipse(const geo_point2& c, const double& r):geo_object(),
		center(c), a(r), b(r)
{
	if(r<=0)
		throw(EME_geometry_logic_error(
				"Ellipse configuration error, r < 0", 6));
	xmin=c.x-r;
	xmax=c.x+r;
	ymin=c.y-r;
	ymax=c.y+r;
	area=M_PI*r*r;
}
geo_ellipse::geo_ellipse(const geo_point2& c, const double& a_in,
		const double& b_in) :
		geo_object(),center(c), a(a_in), b(b_in)
{
	if(a<=0||b<=0)
			throw(EME_geometry_logic_error(
					"Ellipse configuration error, r < 0", 6));
		xmin=c.x-a;
		xmax=c.x+a;
		ymin=c.y-b;
		ymax=c.y+b;
		area=M_PI*a*b;
}
geo_ellipse::geo_ellipse(const geo_ellipse& B) :
		geo_object(B), center(B.center), a(B.a), b(B.b)
{
}
int geo_ellipse::inside(const geo_point2& p){
	double temp=(p.x-center.x)*(p.x-center.x)/a/a+(p.y-center.y)*(p.y-center.y)/b/b;
	if(temp>1){
		return 0;
	}
	else if(temp<1){
		return 1;
	}
	return 2;
}
geo_type geo_ellipse::get_type(){return ellipse;}
}
